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Abstract 

Solutions of Fredholm integral equations of the second kind with oscillatory kernels likely exhibit 
oscillation. Standard numerical methods applied to solving equations of this type have poor numerical 
performance due to the influence of the highly rapid oscillation in the solutions. Understanding of the 
oscillation of the solutions is still inadequate in the literature and thus it requires further investigation. 
For this purpose, we introduce a notion to describe the degree of oscillation of an oscillatory function 
based on the dependence of its norm in a certain function space on the wavenumber. Based on this 
new notion, we construct structured oscillatory spaces with oscillatory structures. The structured 
spaces with a specific oscillatory structure can capture the oscillatory components of the solutions of 
Fredholm integral equations with oscillatory kernels. We then further propose oscillation preserving 
Galerkin methods for solving the equations by incorporating the standard approximation subspace of 
spline functions with a finite number of oscillatory functions which capture the oscillation of the exact 
solutions of the integral equations. We prove that the proposed methods have the optimal convergence 
order uniformly with respect to the wavenumber and they are numerically stable. A numerical example 
is presented to confirm the theoretical estimates. 

Keywords: oscillation preserving; oscillatory integral equation; Galerkin method. 


1 Introduction 

We consider in this paper numerical solutions of Fredholm integral equations of the second kind with 
a highly oscillatory kernel. Specifically, the kernel is a product of a non-oscillatory smooth function 
and a typical known oscillatory function. We assume here that the forcing function is highly oscillatory 
with its oscillation being a combination of harmonic waves. Such an assumption is motivated from the 
plane waves which are used as incident waves in electromagnetic scattering problems. The solutions of 
the integral equations may possess certain oscillation and the error of numerical solutions by standard 
numerical methods may be greatly affected when the oscillation behaves rapidly. Conventional methods 
fail to solve equations of this type. This demands a better understanding of the solutions and devise new 
numerical methods in solving them. 

We choose this type of oscillatory integral equations to study here because they resemble certain 
equations of optics and acoustics [38] and on the other hand they are simple enough to be rigorously 
analyzed to have a good understanding on the oscillation of their solutions. This paper is a beginning effort 
on the understanding of oscillation of the solution of highly oscillatory integral equations and the analysis 
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methods proposed here may supply a potential tool from the mathematical viewpoint to understand highly 
oscillatory phenomena in scientific and engineering applications such as in electromagnetics and acoustics 
| 18 | and in laser theory piiniiii]. 

We remark that several other authors have considered the same integral equation [ini|38]. They 
mainly presented the asymptotic properties of its solution as the wavenumber tends to the infinity. In 
|38j , the asymptotic behavior of the solution was studied with the help of two associated Volterra equations 
and it was proved that the maximum norm of the solution is bounded by the product of the maximum 
norm of the right hand side function of the equation and a constant independent of the oscillation of the 
kernel. A refined result on the asymptotic behavior was presented in m with the assumption that the 
non-oscillatory function in the kernel is smooth and the right hand side function is non-oscillatory and 
smooth. The analysis was conducted through the Neumann series associated with the integral equation. 
Though there is no reference on the numerical solution of this kind of integral equation, two possible 
ways were proposed for the case of the non-oscillatory function in the kernel being 1 in m- One way 
is to solve the integral equations by using the eigenvalues of the integral operator and its corresponding 
eigenfunctions. The shortcoming of this method is that it needs extra cost to compute eigenvalues and 
eigenfunctions. The other way is to approximate directly by the truncated Neumann expansion. This 
method has the poor convergence when the oscillation of the kernel is moderate. Both the methods may 
be difficult to implement, especially when the non-oscillatory factor of the kernel becomes complicated. 

We now review a hybrid numerical method for solving boundary integral equations of the boundary 
value problems of Helmholtz equations (for example, see [12] and references therein). The hybrid method 
combines the conventional piecewise polynomial approximation with high-frequency asymptotics to build 
the basis functions suitable to represent their oscillatory solutions. The idea was hrst proposed in m in 
study of the electromagnetic characterization of wire antennas on or near a three-dimensional metallic 
surface. An ansatz of the high-frequency asymptotics is the base for the hybrid methods and it may 
be derived through the high frequency physical optics (or Kirchhoff) approximation |27] . more precisely 
by the geometrical theory of diffraction [26] or by combining the former approximation with microlocal 
analysis [30] . The ansatz indicates that the solutions of Helmholtz equations may be represented in 
terms of the product of an explicit oscillatory function and a less oscillatory unknown amplitude in each 
zone. The solutions of Helmholtz equations can be obtained by approximating the amplitude with a 
conventional integral method such as the Nystrom method m, the collocation method |23| and the 
Galerkin method nil EH]. Another method in solving Helmholtz equations is the partition of unity in 
which a number of plane waves is introduced on each element in addition to standard piecewise polynomial 
boundary elements [Zl|35|. This method does not require any priori knowledge of the asymptotics of the 
solution which may result in the loss of uniform accuracy with respect to the wavenumber. To maintain 
the effective error, the degrees of freedom need to increase in proportion to the power of the wavenumber 
as the wavenumber tends to the infinity, just as for conventional methods, albeit with a lower constant 
of proportionality. 

There was a recent development in solving the oscillatory Volterra integral equation p Eg im 112]. 
A Filon-type method was proposed in [39] for the numerical solution of the Volterra integral equation of 
the first kind with a highly oscillatory Bessel kernel, based on the fact that its solution has an explicit 
integral expression [39]. A Clenshaw-Curtis-Filon-type method was developed for computing the highly 
oscillatory Bessel transform and was used in solving the oscillatory Volterra integral equation in [42] . 
A Filon-type method and two collocation methods were presented in im for weakly singular Volterra 
integral equations of the second kind with a highly oscillatory Bessel kernel. These methods were designed 
based on the asymptotic analysis of the solutions of the equations. More recently, it was studied in P the 
high-oscillation property of the solutions of the integral equations associated with two classes of Volterra 
integral operators: compact operators with highly oscillatory kernels that are either smooth or weakly 
singular and noncompact cordial Volterra integral operators with highly oscillatory kernels. However, 
numerical analysis of the Volterra integral equations with a highly oscillatory kernel remains a challenging 
problem. 
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Since solutions of weakly singular integral equations has inspired us in developing efficient numerical 
methods for solutions of integral equations with highly oscillatory kernels, it deserves to review some 
related work in the subject. In [8], non-polynomial spline collocations for the Volterra integral equation 
of the second kind with a weakly singular kernel were proposed by making use of the fact that the singu¬ 
larity behavior of its solution can be captured by certain special non-polynomial functions. Singularity 
preserving projection methods were developed in |15j for the Fredholm integral equation of the second 
kind with a weakly singular kernel. These methods gave an optimal order of convergence for the approxi¬ 
mation solutions since the singularity preserving approximation space can preserve the singularity of the 
solutions. Hybrid collocation methods were proposed in m and respectively, for the Volterra and 
Fredholm integral equations with weakly singular kernels. The basic idea of these methods is to enrich 
the basis of a standard approximate space with specific singular functions which can accurately capture 
the singularity characteristics of the solutions. 

In this paper, we make three main contributions to the literature of the solution of oscillatory integral 
equations. We first introduce a new notion to measure the degree of oscillation for an oscillatory function 
and then construct structured spaces with oscillatory structures. The notion of oscillation is defined 
based on the dependence of the norm of a function in a certain space on the wavenumber. It reflects the 
effect of the oscillation of an oscillatory function on the accuracy of its approximation. Structured spaces 
with specific oscillatory structures can capture the oscillation of solutions of oscillatory integral equations. 
The element in such a space may be approximated by an appropriate finite dimensional approximation 
space with its error independent of the oscillation even though the element may be highly oscillatory. 
We then explore the oscillatory property of the solutions of the oscillatory integral equations by using 
the iterated integral operators instead of the asymptotic property of the solutions. The solutions will 
be proved to be in a non-oscillatory structured space with a specific oscillatory structure. Since the 
solutions can be represented by the iterated integral operators, the oscillation of the solutions depends 
on two properties of the iterated integral operators: The non-oscillatory structured space is closed under 
the iterated integral operators and they can reduce the oscillatory degree of oscillatory functions in the 
sense of the new notion of oscillation. Finally, we develop the oscillation preserving Galerkin methods 
(OPGM) to solve the oscillatory integral equations based on the understanding of the oscillation of their 
solutions. The introduction of the methods is inspired by the singularity-preserving methods developed 
for solving the singular Fredholm integral equations of the second |15] and the hybrid methods used in 
solving the scattering problems. The proposed methods have the optimal convergence order uniform with 
respect to the wavenumber and they are numerically stable when the wavenumber is large enough. 

This paper is organized in eight sections. In section 2, we analyze the oscillatory property of the 
solution of the oscillatory integral equation of the second kind and the effect of the possible oscillation 
of the solution on the accuracy of the conventional method in solving the equation. We propose in 
section 3 a notion to characterize the oscillatory degree of an oscillatory function. We then construct the 
corresponding oscillatory spaces and structured spaces with oscillatory structures. We present in section 
4 special results of the structured oscillatory space associated with the Sobolev space. Specifically, we 
show that it is closed under a set of oscillatory Fredholm integral operators. In section 5, we show that the 
solutions of the oscillatory integral equations belong to the structured oscillatory space associated with 
the Sobolev space. In section 6, we develop the OPGM for solving the equation and prove the convergence 
order (uniform with respect to the wavenumber) of the proposed method. We analyze in section 7 the 
stability of a special OPGM which uses the B-spline basis on a uniform partition. In Section 8, we 
discuss the computational implementation of the OPGM. A numerical example is given to illustrate the 
numerical efficiency and accuracy of the proposed method in comparison with the conventional Galerkin 
method. 
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2 Solutions of Oscillatory Fredholm Integral Equations 


We study in this section the oscillatory property of the solution of the oscillatory Fredholm integral 
equation of the second kind. 

We begin with describing the integral equation. Let I := [—1,1]. By C{I) we denote the space of 
continuous complex-valued functions on I and 0(1“^) the space of continuous complex-valued bivariate 
functions on P. Supposing that K £ C{P) and f G C{I) we consider the oscillatory integral equation 


y(s) - ^^y{t)dt = f{s), s £ I, 


( 2 . 1 ) 


where k is the parameter of wavenumber and y denotes the solution to be determined. In highly oscillatory 
problems, k » 1. In this paper, we assume that the range of k is (l,oo). We also assume that K is 
independent of n while / may has some oscillation of wavenumber k. Defining the integral operator JC 

by 

(/Cy)(s) := [ K{s,t)e^^\^-Mt)dt, s £ I, (2.2) 


we may rewrite equation (2.1) in its operator form 

{I-IC)y = f, 


(2.3) 


where X denotes the identity operator. It is well-known that the integral operator K, is compact on C{I) 


and if 1 is not an eigenvalue of JC for any k, then equation (2.3) has a unique solution in C{I) 


We now study how the solution y of equation (2.3) depends on the wavenumber k. It was proved in 


that there exist positive constants c and kq such that for all k > kq 

llylloo < c||/||oo- 

This estimate indicates that the oscillation of the solution y depends completely on the oscillation of the 
right hand side function / in the uniform norm. However, this result does not offer any information about 


the derivatives of the solution y. For the purpose of efficiently solving equation (2.3), we are interested 


in understanding how the derivatives of its solution y depend on the wavenumber k. 

Specifically, we shall bound the Sobolev norm of the solution y above by the power of the wavenumber 
K. To this end, we let P{I) denote the space of functions u on the interval I for which 


\m\2 ■= 

Without ambiguity, we also use the notation 
is. 


\ 1/2 

\u{t)\‘^dt\ < oo. 

1 2 to represent the P norm of bivariate functions, that 


iXllo : — 


r 

I \L{s,t)\‘^dsdt] 


for L £ P{I). 


We next review the notation of the Sobolev space. Let No := N U {0} and for n G N let Z„ := 
{0,I,...,n — I}. For n G No, we let H^{I) denote the Sobolev space {u £ P{I) : pP £ X^(/)} 
with its norm defined by 

1/2 

2 


\u\\h^ := 


E 


,(i) 


'n + 1 


We denote by the space of K-parameterized functions u which satisfy that u £ for any k 

and there exist positive constants c and kq such that for all k > kq, 


\u\\hp < CK^, for all p £ 
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Associated with the operator /C, we introduce two auxiliary operators from L?‘{I) to C{I). Specifically, 
for Lj G C{P),j = 1, 2, and u G we let 

{Ae[Li,L 2 ]u) {s) ■.= J + {—ly J L 2 {s,t)A'^^^~’^^u{t)dt, s€l, (2.4) 

with e G {0,1}. If Li = L 2 , we write L 2 ] as Ae[Li]. When Li = L 2 = K, we have that K. = ^o[Ar]. 

By the Cauchy-Schwarz inequality, we can obtain that 


||»4.e[Li, I/2]u||2 < \J (||Al||| + ||I/2||2) ll'l''l|2, C £ {0)1}- 


(2.5) 


We next present an auxiliary lemma on the derivatives of operators Ae- To this end, for m, n G No we let 
C'[”^’"'](/^) denote the space of functions whose derivatives of order up to m and n with respect to the first 
and second variables, respectively, are continuous. When m = n, we write as CH(/2)^ In the 

next lemma, we let D and use c as a generic constant whose value may change in its appearance. 

Lemma 2.1 Letn G Nq. If K e C^^+I 0 ](^j 2 ) is independent of k andu G ff”(/), then there exist positive 
constants c and kq such that for all k > kq, 


DA A. 




u 


^ CK , j ^ Z„_l_2-i, I G 1in+2: £ — 0, 1. 


( 2 . 6 ) 


Proof: We prove this lemma by induction on 1. The case I = 0 follows directly from (2.5). We assume 
that (2.6) holds foiN < n + 1 and consider the case I + 1. By applying D to Ae[K]u, for e = 0,1, we have 
that 


DiAe[K]u){s) = (A, 


Kim 


u) (z) + iK{A^[K]u){s) + 2eK{s, s)u{s), 


(2.7) 


where e + 1 := [(e + 1) mod 2]. Then by applying to (2.7) with K replaced by with j G 'Ln+i-u 

we find that 


D 


/ + 1 






u) (s) = (.4., 


Kimfi) 


u) (s) + inD^ iA 


e+l 


kHA 


u) (s) + 2eD* 


Ki^’^\g^ s)u(s) 

( 2 . 8 ) 


Using (2.8) with the induction hypothesis, we observe that estimate (2.6) holds for I + 1. Thus, by the 
induction principle, estimate (2.6) holds in general. □ 

Now we are ready to present the bound of the derivatives of the solution by the wavenumber k. 

Theorem 2.2 For a positive integerm, suppose that K gC[™’®](/^) is independent of k and f € 

If y is the solution of equation (2.3), then y G 77^(1). 

Proof: It suffices to prove that there exist positive constants c and kq such that for all p G Z^+i and 

K > Kq, 

WvWhp < ckP. (2.9) 

We prove this by induction on p. When p = 0, we recall an expression in |38j for the solution y 


y{s) = f{s) + J e^ *l{r+(s,t)+r_(s,t)}/(t)dt + J R{s,t] K)y{t)dt, 


( 2 . 10 ) 


where r± independent of k are bounded functions and R{s,t]iPj tends to zero uniformly as k tends to 
the infinity. By applying the if norm to (2.10), we have that 


(i-||7?||2)||y||2<(i + ||r+ + r_||2 


5 



























With the properties of r± and R, there exist two positive constants kq ^ 1 and c such that for all k > kq 

II2/I|2<c||/||2. (2.11) 


The case p = 0 follows directly from (2.11). 


We assume that (2.9) holds for p < m and we consider the casep+1. It suffices to show that there exist 
positive constants c and kq such that for all k > kq, II 2 < ckP+\ Since = {JCy)^P+^H fP+^^ 

and / G it remains to prove that the no rm of {ICy)^P^2} is bounded by ckP^^. Note that 

K G Setting / = p + 1 and j = 0 in (2.6), by Lemma 


2.1 


we have that 


W+^Kyh = \\DP^^ {Ao[K]y) lb < ckP+\ 


Hence, (2.9) holds for the case p + 1. 


□ 


From the inequality (2.11), it is clear that the inverse oiX — K, is bounded by a constant independent 


of K. We state this result as a corollary below. 


Corollary 2.3 If K € C{P) is independent of k and 1 is not an eigenvalue of K, for any n, then the 
inverse of I — 1C exists for any n and there exist positive constants c and kq such that 

sup ||(X-/C)“bl < c. (2.12) 


When the solution of equation (2.3) is highly oscillatory, conventional numerical methods may fail to 
solve the equation. We now elaborate this point by taking the Galerkin methods as an example. Let SJf 
denote the space of spines of order m with h being the maximal distance between two successive knots. 
This space will be defined precisely in Section Let yh G SJf be the approximate solution of equation 
(2.3) obtained by the Galerkin methods. According to [36], if the solution y G then there exists 

a constant /lo > 0 such that for all 0 < h < Lq, 


Wvh - yh < chi' 




By Theorem 2.2, the Ff™ norm of the solution of equation (2.3) may increase in order 0{kI^). As a 
result, there exists a constant c independent of k such that 


hh - vh < c{nhy 


To ensure convergence of the approximate solution y^, we must choose the step size h so that hn < 1. 
When K is large, h will be small. Thus, the resulting linear system will have a large dimension and it is 
computationally costly to solve such a system. This explains why conventional numerical methods may 
fail to the equation when k is large. This motivates us to develop efficient non-conventional numerical 
methods for solving the equation. The difficulty of a conventional numerical method comes from the 
rapid oscillation of the solution. To overcome the difficulty, we are required to understand the oscillatory 
property of the solution. For this purpose, we shall introduce in the next section functional spaces suitable 
for oscillatory functions. 


3 Spaces of Oscillatory Functions 

The main purpose of this section is to introduce a notion which describes the oscillation of an oscillatory 
function. Specifically, we define two kinds of oscillatory spaces, K-oscillatory spaces of order n and 
K-oscillatory structured spaces of order n with oscillatory structures. 

To our best knowledge, there is no appropriate space in the present literature to describe oscillatory 
functions in the context of this paper. Though the space of the functions of bounded mean oscillation 
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(BMO) is relates to oscillation of functions, it can not be used to describe the oscillatory properties of 
oscillatory functions. In fact, it is mainly used to study singular operators. Our first task is to introduce 
appropriate spaces to study oscillatory functions. Then what is an oscillatory function? To have a better 
view on the definition of oscillation, we first review the traditional understanding of oscillation. The 
oscillation was first used in describing the solution behaviour of second order linear differential equations 
|28j and then was extended to describe the solution behaviour of Volterra integral equations. A commonly 
used definition for a solution being oscillatory is given as follows. 

Definition 3.1 A solution x{t) is said to be oscillatory if x{t) has zeros for arbitrary large t (i.e. it has 
infinitely many zeros for t > 0 ); otherwise, a solution x{t) is said to be non-os dilatory 

Though the oscillation defined above is for the solution of some differential and integral equations, it is 
can be extended to the definition of some more general oscillatory functions. 


Definition 3.2 A function y{t) is said to be oscillatory if there exists a known non-oscillatory function 
x{t) such that y{t) — x{t) has zeros for arbitrary large t (i.e. it has infinitely many zeros for t > 0 ); 
otherwise, a solution y{t) is said to be non-oscillatory. 


We note that Definition |3.2| need our intuitive observation of oscillating phenomena to determine a 
non-oscillatory function. According to the extended definition of oscillation, it is defined on an infinite 
domain. Thus it can not give any information of oscillation such as the degree of oscillatory extent on 
finite domain of interest. The influence is that it may be helpless in numerical analysis in the finite 
domain. A simple example is that sinx and sinlO^x on I := [0,1] which are both oscillatory according 
to Definition 3.2 The first function can be approximated numerically using much less computation than 


the second one by standard methods to obtain the same accuracy. Definition 3.2 offers little knowledge 
on this difference. The main problem here is related to another question which is “how oscillatory is an 
oscillatory function?”. 

More recently, high oscillatory problems have been extensively studied by a group of researchers at 
the Isaac Newton Institute of Mathematical Sciences. However, there is no a clearly stated mathematical 
definition on oscillation but it reveals that a function is highly oscillatory if it has rapidly oscillating 
phenomena jUj. It is a conclusion of oscillating phenomena and we can tell which function is more 
oscillating from the oscillating speed, called frequency in time domain or wavenumber in space domain. 
However, the understanding of this oscillation with frequency may still be helpless for numerical analysis 
purpose. It is because in the functional approximation context, what really matter is not just how quickly 
a function oscillates but the effect of oscillation of the function on the accuracy of its approximation. The 
wavenumber alone is not sufficient to describe the effect of oscillation of an oscillatory function that has 
on the approximation accuracy of the function. We now illustrate this point by an example. We consider 
the functions 

2 sin(Kf) 


9j{t) =t + 




t € I, j £ Zt, 


where Zj^ := Tjn \ {0}, for 2 < n G N. Clearly, functions gj, j G Z^ , have the same wavenumber n 
and they oscillate rapidly in the same degree when k is large. However, the effect of oscillation of these 
functions on the accuracy of their approximation are not the same. To see this point, we approximate 
them by using piecewise linear polynomial interpolations on the uniform partition of I with 1281 points. 
We list in Table the maximum errors of the approximation of these functions, where the norm || • ||oo 
is computed by sampling uniformly over I with 2049 points. 

The numerical results indicate that when the wavenumber n doubles, the approximation errors of <71 
increase about four times, those of 52 are about doubled and those of 53 show nearly no change. This 
phenomenon is easily understood when we consider the approximation errors of these functions. For each 
3 e we let gj denote the piecewise interpolation approximation of gj. Then, it can be estimated that 


hj - ffilloo < ch^\\gf\\oo < c/i^(2 + G z|. 


(3.13) 
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Table 1: The maximum errors of the linear interpolations of gj on a uniform partition with 1281 points 


over I. 


K II51-51II00 ||g 2 - g2||oo llffs - gsiloo 


40 

4.89e - 

4 

1.28e- 

5 

9.16e- 

7 

80 

1.95e- 

3 

2.50e- 

5 

9.15e- 

7 

160 

7.80e - 

3 

4.94e - 

5 

9.15e- 

7 

320 

3.11e- 

2 

9.78e - 

5 

9.14e- 

7 

640 

1 .22e- 

1 

1.92e- 

4 

9.09e - 

7 


where c is a constant independent of k and h := 1/1280. It means that for each j G Z^, the error of 


the approximation for gj increases at the speed of as k —)• 00 . Estimate (3.13) indicates that in 

the context of functional approximation, what really matters for the approximation accuracy is how the 
derivative (of the function to be approximated) which bounds the approximation error depends on the 
wavenumber. Specifically, for u G we let 


\u\\c^ := max ||n 
iSZm+l 


(i) I 


and observe that in general the dependence of lluHc"" on the wavenumber k is crucial for the approximation 
of u by piecewise polynomials of order m. 

We now conclude that the traditional understanding of oscillation is a direct description of oscillating 
phenomena and there is no a suitable notion of oscillation which is useful in numerical analysis. This 
motivates us to propose such a notion. The observation made in the above example inspires us to introduce 
the new notion to describe the degree of oscillation of oscillatory functions based on the dependence of 
the norm in a certain function space on the wavenumber. 


Definition 3.3 Let n be a positive number and (X, || • ||x) be a normed space. A function u is called 
K-oscillatory of order n in X if it satisfies 

(1) u is K-oscillatory in X, that is u is n-parameterized and tt G X for any k, 

(2) there exist positive constants c and kq such that for all k > kq 

< c. 

When n = 0, we say that u is non-n-oscillatory in X. 


In the functional approximation context, the space X in Definition |3.3| is normally a Sobolev space 
which appears in the error bound of an approximation of an oscillatory function. Its concrete form 
depends upon the regularity of u, the specihc approximation space and the approximation principle are 
chosen. For example, if the function u is m-times differentiable, the approximation space is chosen to be 
the splines of order m and the approximation principle is the orthogonal projection (resp. interpolation), 
then the appropriate Sobolev space is (resp. IT°°’"*). 

A simple example of K-oscillatory function of order 2 in H^{I) is Based on the definition of k- 

oscillatory, it is clear that the ’’oscillation” appeared in the function is only related to the parameter k. We 
can also understand the K-oscillatory function in some sense by the traditional oscillatory function with 
the wavenumber k. Hence, the parameter k that appears may characterize how rapidly an oscillatory 
function oscillates. In Definition 3.3, the order n reflects the degree of maximum influence that the 


oscillation of a function in space X has on the accuracy of its approximation. The parameter n can also 
be regarded as an index of the amplitude of the oscillatory function in the form of the norm of X. Though 
the concept is proposed to describe the oscillatory functions, some traditional oscillatory functions which 
are independent of k are non-K-oscillatory and some K-oscillatory functions can be non-oscillatory in 
traditional sense according to Definition |3.3[ For example, a non-oscillatory function in traditional sense 







multiplied by a factor k” may be a K-oscillatory function of order n in but not oscillatory according 
to the traditional definition. 

The concept of K-oscillatory functions of order n is a natural development of the classification of 
oscillatory functions in some sence. In the past, we only distinguish the oscillatory functions by different 
wavenumbers or frequencies and now we enrich the classification by distinguishing the difference among 
traditional oscillatory functions of the same wavenumber from the point of their effect on the accuracy 
of the approximation. 

Next, we shall define a K-oscillatory space of order n which gathers all the functions with the same 
maximum influence of the oscillation. 


Definition 3.4 Let n he a positive number and {X, || • ||x) be a normed space. The set '■= {u : u 
is K-oscillatory of order n in X} is called the K-oscillatory space of order n. The space is called a 
non-K-oscillatory space. 


With Definition 


3.4 


and the definition of space in Sectionj^ we have that 

Thus if the hypothesis of Theorem 2.2 are satisfied then the solution y G where X := for 

^ £ T^m+l ■ 

In the following, we shall introduce a kind of important space, K-oscillatory structured space of order 
n with an oscillatory structure. 


|x) be a normed space and M G N. Suppose that 


Definition 3.5 Let n be a positive number, {X, 

:= {sj ■ j G 

the traditional sense. A K-oscillatory structured space of order n with sm is defined by 


is a set of typical known K-oscillatory functions but also oscillatory functions in 


Xfi^n ■ — 



. Uj G G 


(3.14) 


Space Xf^fi is called non-K-oscillatory structured space with the oscillatory structure em- 

It can be verified that for each n, space is a vector space. Every element in the space can be 
represented by the oscillatory structure. For this reason, we call the space the K-oscillatory structured 
space. The choice of the oscillatory structure em depends on the problem under consideration and they 
should describe the main oscillation of interest. For example, when solving equation (2.1), one may 
choose Em as {e*''*, t G /}. 

The non-K-oscillatory structured space with a specific structure is useful in approximating an oscilla¬ 
tory function. In fact, if an oscillatory function is in then we are able to design numerical methods 

to approximate the oscillatory function with the optimal convergence order as we do for a non-oscillatory 
function. This is because we can construct a corresponding structured approximation space with the 
same structure to approximate functions in to avoid the effect from the high oscillation of the os¬ 
cillatory structure. The non-K-oscillatory structured space with a specific structure will serve as a main 


tool in describing the oscillatory properties of the solution of equation (2.1). 


4 Oscillatory Function Spaces Associated with Sobolev Spaces 

We present in this section special results of oscillatory function spaces which defined in Section associ¬ 
ated with Sobolev spaces. These spaces will be used to characterize in the next section the solution of 
the oscillatory Fredholm integral equations of the second kind. 

Through out this paper, let m be a fixed positive integer. According to Definitions |3.4| and 
associated with the Sobolev space H'^{L), for each n G Z^+i, we construct the K-oscillatory space 
of order n and the corresponding K-oscillatory structured space of order n with the 


3.5 
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structure e : t £ I}. In particular, is the corresponding non-K-oscillatory structured 

space. It will be proved that H^q{I) is closed under a set of Fredholm integral operators. 

We begin with an investigation of an oscillatory Volterra integral function related to operator JC. 
We define the oscillatory Volterra integral function below. For a £ I and L £ we define the 

oscillatory Volterra integral function by 

i{s) := f s£l. 

J a 

We recall a well-known result that u £ C^[a, 6], n G N, then 





(-ir 





(4.15) 


where 


n-l 


crn[u]{t) := ^ t £ [a,b]. 


j=0 




(4.16) 


The equation (4.15) can be easily obtained by employing integral by parts. A more general form of (4.15) 
can be found in m- Motivated by (|4.15[), we dehne 


r(s) := -e*'"“crn[L(s, • )](a), s £ I, 


(4.17) 


and 

g{s) := [£(s) — r(s)] s £ I. (4.18) 

Thus, the oscillatory Volterra integral function i has the decomposition 

£{s) = g{s)e'^'^^ + r{s), s £ I. (4.19) 


The expression of r is clear and the property of its derivatives can be easily obtained from those of the 
kernel L. However, the property of the derivatives of g requires further investigation. We next present an 
explicit expression for the derivatives of g. To this end, we first show an auxiliary equality. For j £ Z„_|_i 
and n G No, we let Cn '■= j\(^nLjy, ■ We also need the formula 

[L(s,s)](^’) = s G I, (4.20) 

6=0 

where by [L(s,s)](^) we mean ■^L{s,s). 

Lemma 4.1 If for / G N, L G then for all s £ I 


EEcf(-«) 


l—p^ — lKS 


iN,o)( 


s, s e 


(p-l-q) 


1—1l—l—p 

E E T( 

p =0 5=0 


—IK 


)^-i-P-U(«’P)(s,s). (4.21) 


Proof: Using the Leibniz rule for high order derivatives, we get that 




1 


s e 


p-l-g 


Y. C»_,_^iiKr-^-^- L('?’ 0 )(s,s) 


1 («) 


a=0 


s e 1. 


(4.22) 


We denote by n(s) the left hand side of the formula (4.21). Substituting (4.22) into H and then making 

, sGl. (4.23) 


a change of the summation order of the resulting expression, we have that 
l-ll-l-q 


f ' 





l('?’°)(s,s)' 

[^p=5+a+l 




1 («) 
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By the recursion relation 


Ci = Ci-_\ + Ci_,,j€Zt_,, neN, 


we can rewrite the expression in the curly braces of (4.23) as Hence, (4.23) becomes 

.1 (") 


"W = E e’c, 

q=0 a=0 


, 5 E /. 


(4.24) 


(4.25) 


Equation (4.25) can be expanded further with (4.20). Then we change the summation order of q and of 
a and change q + a to q. This leads to 

1—1 1—1 a. 


= E E E cr,” s), »e /. 

Q;=0 q=OL p=0 


(4.26) 


Changing the summation order again in (4.26) gives that 

i-i q / q \ 


n(s) = EE sei. 

q=0p=0 \a=p / 


(4.27) 


Again, by applying the recursion relation (4.24), we get that Yla=p = Cf ^ and (4.27) becomes 


l-l q 


n(s) = ^^cf“^(-zK)'-^-''L('?-P’P)(s,s), sei. 

q=0 p=0 


(4.28) 


Finally, we change the summation order of p and of q in (4.28) and then in the resulting formula we 
change q — p to q. This yields the desired formula (4.21). □ 

We are now ready to present the explicit expression for derivatives of g. For n,l £ No, let l]n,i '■= 
{{q,p) : q £ I^n,P £ ^i+i,q+p > 1}- 


Lemma 4.2 If L £ C'I"'](/^),n G No, then for I £ Z. 


n+l 


(<?,P)6U„,, P=o 


Proof: By the Leibniz formula we have that 

i 






A direct computation gives 


p=0 


p-i 


s)-Ap\s)\, s£l. 


(4.30) 


£(p)(s) =Y^ \^L^i’^\s,s)e 




q=0 


+ 


J a 


(4.31) 


of r, we obtain that 


By applying (4.15) to the second term on the right hand side of (4.31) and combining definition (4.17) 
itain that 

r 1 i —1i” 

L^'P'^\s,t)e^^^dt-r^P\s)=P'^^(Jn •) (s) + / L^P’^\s,t)e^'^^dt. 


(*'«)’" Ja 


(4.32) 
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To simplify the expression of we substitute the expression (4.16) of an •)] (s) into (4.32) 

and then sub stitu te (4.31) and (4.32) into (4.30). Finally, by combining the simpl ified expression oF^® 
with Lemma 


4.1 


and using the notation i, we obtain the desired formula (^4.29) for 


□ 


We next show that is closed under a Volterra integral operator to be defined below. For a € I 

and L G we define the Volterra integral operator Ca '■ L^{I) —>■ C{I) by 


{Ca[L]u){s) \= / s£l. 


(4.33) 


For each n G No, we define the norm for by ||L||p[n] := maxp^qgz„+i ||oo- Denote by 

non-K-oscillatory space in Thus, for L G there exists a positive constant 

c independent of k such that ||L||g;.[„] < c. 

Lemma 4.3 If L G C^^q{P), then is closed under Ca[L]. 

Proof: In this proof, we shall write Ca[L] as Ca for convenience. Note that for any u G i?™Q(/), there 
exist Vj G G Zj such that 


u{s) = ui(s) + U 2 (s)e*''^ + U 3 (s)e 


s G /. 


Hence, it suffices to prove that for each v G CaV, Ca{ve^^'), and Ca{ve *''') are in 

We first consider the case of a < s. In this case, by the definition of (Cau), we have that 

{Cav){s) = P^^ [\{s,t)v{t)e-^'^^dt. 

J a 

To decompose CqU, we define two functions 

f{s) ■= -e"*'^“crm[T(s, • )v{ • )](a), s G /, 

and 


5(s) : = 


f L{s,t)v{t)e ’'^^dt — r{s) 

J a 


s£l. 


Hence, we have that 

{Cav){s) = g{s) + f(s)e*'^^ s G /. (4.34) 

It remains to prove that g,r & In fact, both g and r satisfy a somewhat stronger condition that 

they belong to Cf^Qp), the non-K-oscillatory space based on C"^{I}. For each I G Zm+i, we apply to 


f and take its maximum norm to obtain that 

m—1 


40 


< 


V — 

^ K<1+ 


, max 

oo ^ s&I 
q=Q 


{LvfP{s,a) 


< 


2'^m 


-max 

K q£Zm 


||l(' 


>'?)! 


oo • 


Since v G L0™g(I), there exists a positive constant c independent of k such that HnUcm-i < c||n||Hm. 
Thus, f G CJ^Qp). By Lemma 4.2, we have that for I G Z^+i 

gii)(s) = - Cf{iky-^-P-pLv)^PP{s,s) + '^Cf{iky-P-^ nLv)^P’P{s,t)e-^^^^-^Ut. (4.35) 


The first summation in the right hand side of (4.35) is bounded by 


(g>p)eu,„,i 


-L-p-q 


max 

se/ 


{LvyPP{s,s] 


< 


2^771^1 


m 


max 


K {q,p)&llm,l 


\\PPP\ 


oo ||n||c'm-l , 


(4.36) 
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and the second summation is bounded by 


p=0 


—p—m 


max 
sei . 


f 

{Lv)^P^'^\s,t) 

1-1 



dt<V22‘^'^+^K^-^ max 

peZi+i 


I G Z, 


m+1 • 


(4.37) 


The combination of (4.36) and (4.37) yields a positive constant c independent of k such that 




<c max ||T^^’'^^||oo||'i’||iy” 

p,(}eZm+i 


We conclude that g G C™q(/). Therefore, CaU G H^q{I). Likewise, we can show that Ca(ne*'''), Ca(ne G 
The case when a > s can be similarly handled. 


□ 


We are now ready to prove the main result of this section. To this end, we define the set A of 
oscillatory Fredholm integral operators by 


A := 


:= |r : r := A[Li, La], e = 0,1, Lj G j = l, 2 } . 

The next theorem reveals that the space is closed under summations of products of operators in 

A. 

Theorem 4.4 If Tj € A, fij G No, az G No and n G N, space is closed under XzLi TYj=i ■ 

Proof: The proof is based on the fact that a Fredholm integral operator can be split into two oscillatory 
Volterra integral operators. Specifically, an operator T G A ca n be split into two Volterra integral 
operators of the same kind as Ca- Hence, according to Lemma 4.3, 77™g(/) is closed under T- 

Suppose that T is a product of any operators in A, that is, T = JlUTj, n G No where Tj G A, j G Z+. 
Since for any u G TjU G 77™q(/), j G Z+, by induction we observe that Tn G 77 ™q(/) which reveals 

that is closed under T. □ 

A direct application of Theorem |4.4| yields the following corollary concerning the iterated operators 
of A. 

Corollary 4.5 If K £ C'H(/2 ) is independent of k, then HfpQ^I) is closed under G N. 

Proof: Since A = Aq[K] and K G is independent of k, the result of this corollary follows directly 

from Theorem 14^ □ 


Corollary |4 .5 1 will play an important role in the next section in understanding the oscillatory structure 
of the solution of equation (2.1). 


5 Oscillation of Solutions of Oscillatory Fredholm Equations 


We study in this section the oscillation property of the solution of equation ( |2.l[ ). This is done by 
considering the iterated integral operators of A. Specifically, the main purpose of this sect ion is to 
establish that if AT G is independent of k and / G the solution of equation (^2.l| belongs 


The solution y of (|2.1|) can be represented by using iterated integral operators of A. In fact, by 

(5.38) 


successive substitutions, we obtain from (2.1) that 

n—1 


y = A"y + j;AV, n G N. 
i=o 
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It is clear that the oscillatory property of the solution depends on the properties of the iterated integral 

that if / G Yl'jZo f ™ remains to 

the solution 


4.5 


operators of /C. It follows from Corollary 
prove that there exists a number re G N such that fC^y G According to Theorem 


2.2 


of (2.1) belongs to Hence, it suffices to show that there exists an re G N such 


that /C- : ^ 

The next lemma shows that for L G satisfying an additional condition, the integral operator 


Ca[L] defined by (4.33) can ease oscillation of functions in 


Lemma 5.1 Let re G ^^+ 1 - V L ^ and there exists a positive constant cq independent of k 

such that 

0<a + /3 + l<re, (5.39) 


max 

se/ 


f/iere C,[L] : F-(/) ^ 

Proof: We provide proof for the case a < s only since the other case can be similarly handled. According to 


the definition (4.33), we have that {Ca[L\u) (s) = rei(s)e*''^, for s € I, where rei(s) := L{s, f)e ^'^^u{t)dt. 


It suffices to prove that rei G 
Differentiation of ui leads to 


u 


I’w-E E'cf.j.jLi".”)!*,, 

p=0 5=0 


1 


/ L^^'^\s,t)e-^'^^u{t)dt, f G Z^+i. 


(5.40) 


By using the Cauchy-Schwarz inequality, the norm of the second term on the right hand side of (5.40) 
is bounded by 2 maXgg 2 ™ > i 


'm+1 


Pl|2- 


We next bound the first term (denoted by di{z)) of the right hand side of (5.40). Since u G 


by the Leibniz formula, it is easy to verify that there exists a positive constant c independent of re such 
that 

r .. 1 ini 

< cnf, P G 

2 



(p) 




(5.41) 


The L?' norm of the terms with p <l — nui5iis bounded by ck^ ^ according to (5.41) and L G cj^ If 

p > l — n, then the factor (s, s)] ^ be represented by a linear combination of (s, s), 

0<Q!-|-/3<f — 1— p. Thus, the norm of the terms in Si{z) for p > I — n is bounded by 


according to the condition (5.39) and (5.41). Hence, there exists a positive constant c ind ependent k 
I 2 < for all I G Zm-i-i- Combining the bounds of the two terms in (5.40), we obtain 

that there exists a positive constant c independent re such that ||re ^^^||2 < cmax{K^“"', 1}, for all I G Zm-i-i- 

□ 


such that 
that there 
Therefore, we have that rei G 


The oscillation order of a function in H^{I) is m and that of a function in is m—n according 



Lemma 5.2 Let re G IfLj , j = 1, 2, are the less oscillatory kernels of order re, then for e = 0,1, 

A[Zi,L2]://ru)^^K”:m-nW- 


Proof: Note that the operator Ae[Li, L 2 ] can be split as two Volterra integral operators of the type Ca[L]. 
The result of this lemma follows directly from Lemma 5.1 □ 
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We next study the effect of the iterated integral operators (HO) of JC applied to This is done 

by using induction on the order of the iterated integral operators and finally choosing an suitable n G N 
such that K 7 ^H'^{I) C Notice that the iterated integral operator of K. of order 1 is equal to fC. 

Since K G is independent of k, we conclude that K G C^^q{P) and there exists a constant c 

indepen den t of k such that max^g/ \K{z,z)\ < c. Hence, K is the less oscillatory kernel of order 1. By 
Lemma 5.2, JC maps into This means that JC reduces the oscillatory order of functions 

in F-()T^y 1 . 

We shall show that for n G reduce the oscillatory order of functions in by n. 

Because the IIO of a high order is the composition of K, and the IIO of one order lower, the property 
of the IIO of a high order can be derived through induction on the order. To decompose a composite 
operator, we need to decompose its corresponding kernel. Since the kernel of the composite operator is 
represented by an oscillatory integral whose integrand has the non-differentiable oscillator, we need to 
split the integral domain. To this end, we introduce auxiliary functions to represent the decomposition 
of the kernel. For e = 1, 2, we let := K. For 1 < n G N and e = 1, 2, we define 


Ln-l,e{s,t,x) := K{s,x)Kn-l,eix,t), {s,t,x) G 


(5.42) 


The function Ln-i,e is the less oscillatory part in the integrand of the kernel of the composite operator. 
To represent the decomposition, for 1 < n G N and e = 1,2, we define 


Jn,e{s,t) := (-1) 


e-l 


Ln-l,e{s,t,x)dx, (s,t)G/^. 


(5.43) 


Let V := {(l,-l,2,t), (2,-l,2,s), (3,1, l,s), (4, l,l,t)}. For 1 < n G N and {j,j,ej,z) G V, 

Qn,j{s,t) :=-j [ Ln-i,ej{s,t,x)e‘^"^'^'^dx, {s,t)£l^. (5.44) 


Noting that Qn,j still possesses certain oscillation, with the decomposition (4.17) and (4.18) of oscillatory 
Volterra integrals, for 1 < n G N and {j,j, ej, z) G V, we define 


:—je dm [L^— (s, t, •)] (j), (s, t) G I , 


Gn,j{s,t) := {Qn,j{s,t) - Rn,jis,t))e 


— 2ijKZ 


{s,t) G 


(5.45) 

(5.46) 


where dm is dehned in the same way as am with k being replaced by 2jK. It is easy to check that 
Rn,i = Rn,2 and = Rn,A and for {j,j, ej, z) G V, 


Qn,j{s,t) — Gnj(s, + Rn,j{s,t). 

For 1 < n G N and e = 1, 2, we let 

Rn,e{s, t) := Gn,e{s, t) + Jn,e{s, t) + Gn,e+2{s, t), (s, t) G 


(5.47) 


(5.48) 


The functions Kn,e-, e = 1, 2 are the less oscillatory part of the kernel for the new integral operator which 
is separated from the composite operator. To better understand the functions defined here, readers are 
referred to the proof of Lemma 5.6 We shall present the property of Kn^e, n G '^m+ii e = 1) 2. To this 
end, we study in the next lemma the property of Jn,€-, for e = 1, 2 and 1 < n < m. 

Lemma 5.3 Let 1 < n < m. If K e C'H(/ 2 ) is independent of k and ‘I o-xe the less 

oscillatory kernels of order n — 1, then Jn,e, e = 1 , 2 , are the less oscillatory kernels of order n. 
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Proof: Since K G is independent of k and Kn-i,e £ = 1,2, it is straightforward to 

show that Jn,e £ e = 1, 2. It suffices to prove that Jn,e, e = 1, 2, satisfy (5.39). To this end, we 

consider the maximum bound of Jn’e\s, s) ioi I + p < n — 1. For the concise representation of Jnf\s, s), 
let Sa{s,t) := (s,t, s), a G and Pp{s,t) := (3 G Zp. A direct calculation with the 

help of (4.31|) yields 


i-i 


p-i 


4f{s, s) = (-1)'-' ^ o) - E P/ 


(0,p-l-/3) 

/3 


(s,s) 


(5.49) 


L a=0 


3=0 


By using the Leibniz rule, we have that 

l—l—a 


7=0 




and 


7=0 


P-1-/3 


(p- 1 -^- 7 ) 


(5.50) 


(5.51) 


t=S 


Since K G CH(/7 is independent of k, the derivatives of K appeared in ( |5.50[ ) and ( |5.5l[ ) are uniformly 
bounded by a constant c independent of k. Note that the orders of derivatives of Kn-i,e in (5.50) and 
(5.51) are at most n — 2. Recalling that K^-ieP = 1)2 are the less oscillatory kernel of order n — 1, 


the derivative of Kn-i,e in (5.50 


^_ and (5.51) is uniformly bounded by Thus, Sa ^ °‘’^\s,s) 

and Pg^’^ ^ s) are bounded uniformly by for a positive constant c independent of k. By 


□ 


using (5.49) we conclude that Jn,eP = 1)2, satisfy (5.39). 

Next, we consider functions Rnj and Gn,j- 

Lemma 5.4 Let 1 < n < m. If K G C'H(/7 is independent of k and Kn-i,eP = 1)2, are the less 
oscillatory kernels of order n— 1, then for j G Z^, Rnj,Gn,j G and Gn,j are the less oscillatory 

kernels of order n. 

Proof: Since K G ^^(/S) is independent of k and Kn-i,e £ e = 1, 2, by a proof similar to that 

of Lemma 4.3, Gn,j and Rnj, j G Z^ can be proved in 

We next show that Gnj, j G Z^ are the less oscillatory kernels of order n. Because Gnj, j G Z^ have 
a similar structure, they can be proved in the same way. Thus, we verify only the case of Gn,i- Suppose 
that l,p G N and ( +p < n — 1. From (5.46), Gn,i has an explicit expression 


Gn,l{Sjt) — 


Ln-l,2{s, t, x)e + e^^'^dm [Ln-l, 2 (s, L •)] (-1) 


1-1 


e^^^\ {s,t)Gl\ 


For convenient presentation, we let Ss{t,x) := L^’_^fl{s,t,x). Then, the Lth order partial derivative of 
Gn,i with respect to the variable s has the form 




SsiL x)e-2-"dx + e^^^dm [Ss{t, • )] (-1) 


2iK, ; 


f-1 


^2iKt 


(5.52) 


Applying Lemma 
'lilp) 


4.2 


Up) 


to (5.52) leads to the following explicit expression of 

p 


Cll'P>{s,t) = - G^{2iK)P-^-°-^Si>^’^'>{t,t) + '^Cl{2iK)P-^-"^ f 

{a,l3)eVm,p h=0 

(5.53) 


16 




























To show that Gn,i is the less oscillatory kernel of order n, we consider the bound of We 

consider the two summations on the right hand side of (5.53) separately. We first consider the first sum. 
The terms satisfying a + j5 > n — 2 m. the sum are bounded by since p—1 — a —/3 < p+1 —n. For 

any a + /3 < n — 2 , (s, s) is bounded by due to the fact that Kn- 1^2 is the less oscillatory 

kernel of order n — 1. Thus (s, s) is uniformly bounded by ” for all a, /3 G Um,p. 

Therefore, the first sum on the right hand side of (5.53) is bounded by 

We next consider the second sum. To this end, we let U := Um,p U {(a,/3) : a = m,/? G "^p+i} and 
observe that 


(a,^) G U. 

7=0 


(5.54) 


Since K G is independent of n and Kn- 1,2 £ there exists a constant c independent of k 


such that ||iF(^’'>')||oo < c and \\K^ 


(«-7.^) I 
n-1,2 I 


< c for alH < m, a < m, /3 < m and 7 < a. These estimates 


together with equation (5.54) lead to the estimate ||oo < 2“c^. It follows for n > 1/2 that 


'£cl{2tK 


,p-( 3 -m 


/3=0 




< ^C^(2n)P-^-™2™+ic2 < 2^p+^c^kP-^. 

/3=0 


Thus, we conclude that there exists a constant ci independent of k such that the second term on the 
right hand side of (5.53) is bounded by 

We now combine the two cases above to obtain the desired result of this lemma. Noting that p +1 — n 
and p — m are both no more than I +p+ 1 — n, we have that (s, s) is uniformly bounded by 
with respect to s for a positive constant c independent of k. Consequently, Gn,i is the less oscillatory 
kernel of order n. □ 

In the next lemma we establish that Kn^e is the less oscillatory kernel of order n. 

Lemma 5.5 If K £ is independent of k, then Kn,e, e = 1,2 are the less oseillatory kernels of 

order n, n & ^m+i- 


Proof: According to the discussion presented after Lemma 5.2 AT is the less oscillatory kernel of order 1. 
The result of this lemma may be proved by induction on n with a help of Lemmas |5.3| and 5.4 □ 


We need a set of auxiliary integral operators to represent the decomposition of the IIO. For n G N, 
we define /C„ := .4.o[A'„,i, Kn, 2 \- In particular, we have that /Ci = /C. We present in the next lemma the 
decomposition for iterated integral operators of high orders. 


Lemma 5.6 If K £ is independent of k, then for u G HJf{I) 

KPu = JCnU + u, u G H^q{I), n G 


(5.55) 


Proof: The proof of this lemma is done by using induction on re. The case for re = 1 trivially holds. 
Suppose that the case for re < rre is proved and we consider the case re + 1. By applying JC to equation 
(5.55), we have that 

/C^+^re = K {Knu + u) = ICICnU + ICv. 

Since v G H'^q{I), it follows from Corollary 4.5 that ICv G HI^q^I). It remains to consider the composite 
operator JCJCn- We decompose its kernel which has the form 


K{s, x)Kn{x, 


/-I 
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where Kn{x,t) := Kn,i{x,t) if x > i and Kn{x,t) := Kn^ 2 {x,t), otherwise. We next decompose the kernel 
by splitting the integral domain to remove the absolute signs. We consider two cases t < s and t > s. 
When t < s, the integral domain is split into three parts [—l,t], [t, s] and [s, 1]. In this case, the kernel 
equals 

t) + t) + e-*"(^+*)g„+i,3(s, t). (5.56) 

By substituting (5.47) into ( 5.56| ), with the definition of we write the kernel as 

+ + 72^+1,3(5, t). 

In the other case, the integral domain is split into [—l,s],[s,t] and [t, 1] and accordingly the kernel 
becomes 

We have already used that Rn+ 1,2 = Rn+ 1,1 and Rn+i,i = 7?n+i,3- By introducing the function 


ui(s) := J i?n+iq(s,t)e* + J Rn+i,3{s,t)e u{t)dte~ 

according to the decomposition of the kernel, we have that 

^^ 72 ^ - + ^ 1 * 

Since for j £ Z^, Rn+ij £ as proved in Lemma |5.4| 


HJ^q{I). Hence the decomposition (5.55) for is verifiec . 


it is straightforward to show that vi £ 

□ 


We are now ready to show a property of the IIO of /C, which is crucial for establishing the main result 
of this section. 

Theorem 5.7 If K £ is independent of k, then /C" : 77™(/) —)• n £ ’^^+ 1 - 


Proof: Since K £ is independent of k, by Lemma |5.6l we have the decomposition that IC^u = 

JCnU + V, V £ for u £ According to Lemmas 5.5 and 5.2, we observe that JCn ■ —)• 


Noting that C we conclude that K'^u £ for u £ 


□ 


Theorem 5.7 implies that KP can reduce the oscillation order of oscillatory functions in In 

particular, when n = m, we have that /C™'77™(/) C 77™g(/). This is our desired result. As a direct 
consequence of Theorem 5.7 we obtain the main result of this section, which concerns the oscillatory 
property of the solution of equation (2.1). 


Theorem 5.8 Let y be the solution of equation (2.1). If K £ is independent of k and f £ 

then y £ (I). 


Proof: According to Theorem 2.2, we hav e that y £ HJf jl). Since K £ is independent of k and 

/ G HILq{I), by employing Corollary 


from equation (5.38) with n = m. 


4.5 


and Theorem 


5.7 


the result of this theorem follows naturally 

□ 


Theorem 5.8 characterizes the oscillatory structure of the solution and it serves as a base for developing 
efficient numerical methods for solving the equation (2.1). 
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6 Oscillation preserving Galerkin method 


In this section we shall develop an oscillation preserving Galerkin method (OPGM) for the numerical 
solution of equation (2.1). Specifically, we shall introduce to the usual spline space certain simple, basic 
oscillation functions which characterize the oscillation of the solution of equation (2.1). We shall establish 
the convergence theorem for the proposed OPGM. 

We begin with describing the usual spline space which is used to approximate the solution of the 
equation. Given a partition of /, A : —1 = to < h < ■ ■ • < < ^N+i = 1; let h := max^.gg+ {tj — tj-i) 

and Ij := {tj-i,tj), j G ^iv+ 2 - denote by 11^ the set of polynomials of degree less than m. Let 
5™’^(A) := {re G : w\i- G 11^,^ G Z^_|_ 2 } where (when ly = —1, the elements of 

this space are allowed to possess jump discontinuous at their knots). This space is called the space of 
spline functions of degree m — 1 with knots at ti,t 2 , ■ ■ ■ ,t]s[ of multiplicity m — 1 — v and its dimension 
is given by d := m{N + 1) — A(1 + v). The smoothest space of nondegenerate splines is the one with 
v = m — 2 which has dimension N m. Let S'™ := S™’™~^(A). The space S™ has the B-spline basis 
{Bj : j G Zrf} of order m where d = N + m. Let denote the orthogonal projection from onto S™. 
It is well-known that V'j^ has the properties: 

(i) For all u G — "^112 —>• 0 as h —)■ 0. 

(ii) There exist constants c > 0 and ho > 0 such that for all 0 < h < Lq, W'P'f^u — n ||2 < ch'^\\u\\H^, 
for all u G iL™(I). 

As indicated at the end of Section]^ approximate solutions of equation (2.1) from standard spline spaces 
suffer from the oscillation of its exact solution. Theorem 5.8 suggests that we should include in the usual 
spline space additional simple, basic oscillation functions which characterize the oscillation of the solution 
of the equation. This gives the OPGM which we describe next. 

According to Theorem 


structure {e* 


5.8 


we shall use an extended spline space based on space S'™ with the 


: s € 1} as our approximate space for the numerical solutions. Namely, we let 


■— -|-^26 + Wo : Wj G G Z^} (6.57) 

and call the structure space. Let denote the orthogonal projection from onto Since 

this projection maps a function in into with an oscillatory structure, we call it an oscillation 
preserving projection, and call the solution yh of the equation 


Vh = VhlCyh + Vhf, 


(6.58) 


or equivalently 


{Vh, Vh) = {JCyh, Vh) + (/, Vh), for all Vh G V^h, 


(6.59) 


an OPGM approximate solution in for equation (2.1). Since the solution y of equation (2.1) and 
share the same oscillation structure, yh can capture the main oscillation of the solution y. We expect 
that the numerical solution yh approximates the y in an optimal order. 

In order to analyze Vh, we define two K-oscillatory spaces IP™ := {we'‘ 

Gorrespondingly, we define two linear operators from onto kP™ by 

2 


: w G 5™}, for r = —1,1. 


{V'[u){s) := {Vhw){s)e'^'^'^^, s e I where w{s) := u{s)e for all u ^ Lr 


(6.60) 


Note that {{Vl)‘^u){s) = (P,"(P»e"-)(s) = {V'^{V'^w)){s)e^^-^ = (P»(s)e"- = {Vlu){s),s G /, for 
all u G L‘^{I)- Operators V]) are projection operators from onto kP™, r = ±1. 

We recall one known result taken from [T]. 

Proposition 6.1 Let X be a Banach space. If T,Tn '■ X ^ X are bounded linear operators with Tn ^ T 
pointwise as n —>■ oo, then 

\\{Tn - T)C\\ 0,n ->CX) 

for each compact operator C : X ^ X. 
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To show the error bound of OPGM, we also define for u G 

1/2 


U 


HU 


:= inf < 



: u = wie"^'+W2e + W3,Wj e H'^q, j = 1,2,3 


We are now ready to prove the main theorem of this section which gives the convergence of OPGM. 

Theorem 6.2 Suppose that K G is independent of k, f £ and 1 is not the eigenvalue 

of JC for any k. If y is the solution of equation (2.1), then for sufficiently large k, there exist positive 
constants c independent of n and ho > 0 sueh that for all 0 < h < ho, the OPGM has a unique solution 
yh £ ^ffh satisfies 


\\y - Vhh < ch^WvW 


H. 


K,0 


(6.61) 


Proof: From (2.3) and (6.58), we obtain that 

(X - P/j/C) {yh -y)= VhV - y- (6.62) 

Since Vh is the orthogonal projection from onto we have that \\Vh\\ = 1- By property (i) of 
and the inequality 

\\PhU - U\\2 < (1 + W'PhlDW'P'hU - U\\2, 

we observe for u G L‘^{I) that \\VhU — u ||2 —)• 0, as /i —)■ 0. Hence, by employing Proposition 6.1 together 
with the compactness of /C, we have that ||Pfe/C — /C||—>-0, as/i—)-0. The existence of {X — K.)~^ ensures 
that there exists an /iq > 0 such that for all 0 < h < ho, 

1 


Noting that 


we obtain by inequality (6.63) that 

\\{X-VhlC)-^\\ < 


\\{z- ic)-\Vhic -m < 

X - VhIC = {X-1C)[X+{X- /C)-1(/C - VhlC)] 

ii(x-/c)-i 


(6.63) 


<2||(X-/C) 


-ii 


l-||(X-/C)-i(Ph/C-/C)|| 

This ensures the existence of (I — VhK,)~^■ Thus, we have from (6.62) that 

y-yh = {X- VhfC)~^{y - VhV)- 
It remains to estimate \\y — yh\\ 2 - Theorem 


(6.64) 


(6.65) 


5.8 


ensures that y G H(P‘q{I). This means that for 

( 6 . 66 ) 


sufficiently large k there exist Wj G H)f'Q{I),j G such that 

y{z) = t(;i(s)e*''® + W 2 {s)e~'"^^ + W 3 {s), s £ I. 

Combining ( |6.64 ), ( |6.65 ) and ( 6 . 66 ) together, we derive that 

Wy-Vhh < 2 ||(X-/C)“^||(||u;ie*'^' - P/i(rcie *'"')||2 + ||u; 2 e“*'"' - P/i(tC 2 e “*'"')||2 + ||tC 3 - P/itcsIb)- (6.67) 

According to Corollary |2.3t for sufficiently large k there exists a positive constant c independent of k such 
that ||(X — /C)“^|| < c. It suffices to consider the three terms in the parenthesis in (6.67). By definition 
de:^ , we obtain that 

\\wje"'^' - Pft(tcie *^^''"')||2 < 2\\wje’-'"' - Vl^ {wje "'^^'^')\\2 = 2\\wj -'P'hWj\\ 2 , ( 6 . 68 ) 

where tj = j — 2, j £ Z^ and P° := P)^. The desired result (6.61) is obtained by substituting ( 6 . 68 ) into 
(6.67) and using property (ii) of P)^. 


□ 


Noting that for y £ HffQ{I), 




K,0 


is bounded by a positive constant independent of n, solving the 


oscillatory equation (2.1) by OPGM gives the optimal order of convergence which equals the order of the 
spline functions used in the approximation. 
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7 Stability Analysis 


In the section, we investigate the stability of OPGM based on the B-splines, in terms of the condition 
number of the coefficient matrix of the linear system that results from the method. We shall establish 
that when k is large enough, the condition number of the discrete system is bounded. 

We first describe the linear system for OPGM based on the B-splines. Suppose that space SJ^ is 
defined based on a uniform partition and the B-spline basis {Bj : j G Z^} of order m is chosen for space 
5™. By the definition (6.57) of we have that B := : j G Z^} is a basis for the 

structure space bTh. Tne OPGM approximate solution of the equation can then be represented by 


yh{z) = '^aijBj{s)e 
i=i 


d 

^^^ + T.a2,jB, 


i=i 


s)e 


j=i 


s G /, 


where coefficient vector a := [aij : I G Z^,j G Z^]^, is determined by the linear system obtained by 


substituting the above expression of yh in (6.59) and replacing Vh by each basis function in B. To set 
up the linear system, we let Cp ■.= p — 1, for p G Z3, and for p G Z3, we dehne the vectors fp := 
[(/, : j G Zrf] and use them as blocks to form the vector f := [fp : p G Z3]. For p, (7 G Z3, we 

introduce matrices 

Ep,g := [{Bie^^^^-,Bje^^^^-) : l,j G Z^] , Kp,, := [{JC : l,j G Z^] 


and define := [Ep^g ■ P,(l ^ ^ 3 ] and := [Kp^g ■ P,Q ^ ^ 3 ] • Thus, the discrete linear system (6.59) of 
OPGM can be rewritten in the matrix form 


{Ed - Krf)a = frf. 

We next study the condition number of matrix (E^ — K(i). For this purpose, we recall two well-known 
facts. 

The first fact concerns the condition number of the coefficient matrix of a general Galerkin method, 
which was presented in [ 2 ] (pp. 94-97). Suppose that the operator equation (A —/C)y = / is solved by the 
Galerkin method with the basis {(j)Q ,..., (l)d-i} using the orthogonal projection operator Vn- The resulting 
approximate operator equation is given by {\ — 'PnK,)y = Vnf ■ We use to denote its coefficient matrix. 
Let Fn := \{(pp(t>j) : l,j G denote the Gram matrix of the basis. Then the condition number of the 
coefficient matrix A„ satishes the bound 


Gond(An) < Cond(r„)Gond(A — PnAf), 


(7.69) 


where Gond(A) denotes the condition number of the matrix A. 

The second fact regards the condition number of the B-spline basis of order k HE]. We denote by 
t := {tj : j G "^d+k} the equal-space sequence of knots with the first k and last k knots in t being the 
same and to = —1, td+k-i = 1- Let J\f := {Nj : j G Z^} denote the B-spline basis of order k on t, that is 

Nj{t) := {tj+k - tj)[tj,.. .,tj+k]{ • - t)+“\ j ^^d,t£ I, (7.70) 


where [tj ,..., tj^k]u denotes the k-th order divided difference of u at the nodes tj,..., tj+k and a+ := 
max{Q!, 0}. The condition number of the basis {A^} is defined as 


Gondfc := sup 

t 


Ib| 


sup ■ 
b^O 


E 


IL 

bjNj 


sup ■ 
1 2 b^O 


EjeZd bjNj I 
l|b||;2 


It is proved that there exists a positive constant Dk which only depends on k such that. 


Gondfc < Dk- (7.71) 

We denote by G := [{Ni,Nj) : I, j G Z^] the Gram matrix of the B-spline basis. The following lemma 
shows that the condition number of G is bounded by a constant. 
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Lemma 7.1 If Af is the B-spline basis of order k defined by (7.70), then Cond{G) is bounded by a 
constant independent of the order of the matrix. 


Proof: We shall connect Cond (G) with Cond^ and make use of estimate (7.71). Since matrix G is 
symmetric, we have that 

Cond(G) = An,ax(G)/A„,in(G). (7.72) 

We next associate Cond(G) with the condition number of the basis Af through the maximum and min¬ 
imum eigenvalues. For this purpose, we derive an alternative expression for the Cond^. Recall that 




= b^Gb. It is straightforward to have that 


b^ Gb 

sup ||, 112 
b^O l|t>ll;2 


sup Xj 


b^O 




E 




— Aniax(G) 


and 


||t>||p 1 / ■ r b^Gb 1 / • r \ 

b^ob^Gb 'b^ollbllf^ 


l£l.d V 


— l/Amin(G). 


These formulas together with the definition of Cond^, we get that 

Condfc = sup \/Amax(G)/Amin(G). 


It follows from (7.72) and (7.73) that 


Cond(G) < Cond|. 


(7.73) 

(7.74) 
□ 


We thus obtain the desired result from (7.71) and (7.74). 

We now present the theorem regarding the condition number of the matrix — K^. 

Theorem 7.2 Let K G C'[™'1(/^) he independent of n, f € and y be the solution of equation 

(2.1). If the equation \2. j| ) is solved numerically by OPGM based on the B-spline basis with the uniform 
partition, then there exist constants c and ho such that for all 0 < h < ho there exists a constant kq such 
that for all k > kq, Cond{Ed — K^) < c. 

Proof: According to the bound ( |7.69 ) of the condition number of discrete system of general Galerkin 
methods, we have that 

Cond(Ed - Kd) < Cond(Erf)Cond(X - VhlC). 


By Corollary 2.3 there exist constants ci and ki such that for all n > ki, Cond(X — 1C) < ci. Since 
I — VhIC —)• X — /C, as /i — 7 - 0, there exists a constant ho such that for all 0 < /i < /iq and k > ki, 
Cond(X - Rfe/C) < Cl. 

Next, we bound the condition number of E^ for all h < ho. We write Eq^o as E and observe that 
Ep^p = E, for p G Z 3 . When p,q G Z 3 , p ^ q, hy the Riemann-Lebesgue Lemma, the element in Ep^g 
tends to 0 as K —)■ 00 and thus, we get that 


lim Eld = Diag(E, E, E) 

k;^oo 


(7.75) 


We notice that E is the Gram matrix of the B-spline basis of order m which is independent of k. Lemma 

(7.76) 


7.1 ensures that there exists a positive constant C 2 independent of h and k such that 

Cond(E) < C 2 . 
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Due to the fact that Diag(E, E, E) and E are positive-definite matrices, we obtain that 


Cond(Diag(E,E,E)) = 


Amax (Dia-g(E, E, E)) _ Amax(E) 
Amin (Diag(E,E,E)) 


= Cond(E). 

Amin 1 


(7.77) 


Combining (7.75), (7.76) and (7.77), we have that there exists a constant Kh such that for all k > 


Cond(Erf) < C 2 . We note here that the value of depends on h. By setting c := ciC 2 and kq := 
max{Ki, K/i}, for all k > kq, Cond(Erf — K^) is bounded by c. □ 


Theorem |7.2| ensures that OPGM based on the B-spline basis is numerically stable when the wavenum¬ 
ber is large enough. 


8 Numerical Experiments 

In this section, we present a numerical example to demonstrate the approximation accuracy, the order 
convergence and the stability of the proposed method. We also compare it with the conventional Galerkin 
method (CGM). 

In this example, for simplicity we choose K(s,t) = 1. For the numerical comparison purpose, we 
choose / as 

f{s) = (-]s^ + (1 - —)s^ + + -^s - ^ + 


r,2in 


+ Qik + — 4zk^) — -e*'^ 

8k^ k 


. 2i 

-|- 1-, s G I 

K 


SO that the exact solution of equation (2.1) is given by y{s) := 1 -|- , s G I. It is straightforward to 


obtain the norm of y which is ||y ||2 = for any k. 

We choose the uniform partition A := {tj := —1 -|- jh : j G ZAr_|_ 2 } for interval I, where h := 
with t-i = to and tiv +2 = t 7 V-i-i- We choose m = 2 and the B-spline basis {Bj ^2 ■ j £ '^N+ 2 } of order 2 as a 
basis for S'^. That is, Bj ^2 is a piecewise linear function whose support is and Bj^ 2 {tj) = 1, for 

j G ZAr_|_ 2 . The CGM will use the approximation space S'^ while the OPGM will use the non-K-oscillatory 
structured space based on S'^. For a fair comparison, we shall choose different N values for the two 
methods. We use N to represent the order of the coefficient matrices. We have that N = 3{N -|- 2) for 
the OPGM and N = N + 2 for the GGM. To make the dimension of discrete matrices comparable, the 


values of N for the GGM are set four times of those of N for the OPGM. From Theorem 6.2, we have that 
\\y — yh \\2 = 0{h?‘) where y is the exact solution of equation (2.1) and yh is the OPGM approximation in 


The relative error in the L? norm is evaluated numerically by the trapezoidal quadrature rule with 
2048 uniformly distributed notes, which is 


\y - yh \\2 


l|y|| 


e-N 




1/2 


\ ^ 2048 
\i=i 


\y{sj) - yhisj)^^ 


(8.78) 


where Sj = — 1 - 1 - 
defined by 


for j G Z+ 


1024 


2049- 


We use C.O. to denote the computed convergence order, which is 


G.O. := log f—) /log(2). 

\e2N J 


Likewise, we let denote the relative error of the approximate solution generated by the GGM and its 
computed convergence order is defined accordingly. 


We perform two experiments. In experiment one, we solve equation (2.1) with fixing the wavenumber 


K to 50,500, 5000 and 50000 and varying the order N of the coefficient matrices. We report the numerical 
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Table 2 : The related errors of the CGM (Left) and the OPGM (Right) with different values of k 


(a) K = 50 (b) K = 500 


N 

Bn 

C.O. 

N 

Cat 


C.O. 

N 

Cat 


C.O. 

N 

Bn 

C.O. 

66 

4.23e- 2 


54 

1.36e- 

3 


66 

2.50e - 

1 


54 

1.67e- 3 


130 

1.08e- 2 

1.98 

102 

7.24e - 

4 

0.92 

130 

2.49e - 

1 

0.00 

102 

4.54e-4 

1.88 

258 

2.75e- 3 

1.97 

198 

1.46e- 

4 

2.31 

258 

2.38e - 

1 

0.06 

198 

1.20e-4 

1.92 

514 

7.03e - 4 

1.97 

390 

1.25e- 

5 

3.54 

514 

7.60e - 

2 

1.65 

390 

2.67e - 5 

2.16 

1026 

1.93e-4 

1.87 

774 

1.20e- 

6 

3.38 

1026 

1.93e- 

2 

1.98 

774 

1.98e- 5 

0.43 

2050 

5.20e- 5 

1.89 

1542 

1.17e- 

6 

0.04 

2050 

5.18e- 

3 

1.90 

1542 

3.06e - 6 

2.69 



(C) K = 

5000 







(d) K = 


50000 



N 

Cat 

C.O. 

N 

Cat 


C.O. 

N 

Cat 


C.O. 

N 

Bn 

C.O. 

66 

2.50e- 1 


54 

1.67e- 

3 


66 

2.50e - 

1 


54 

1.67e- 3 


130 

2.50e- 1 

-0.00 

102 

4.55e- 

4 

1.87 

130 

2.50e - 

1 

-0.00 

102 

4.55e-4 

1.87 

258 

2.50e- 1 

-0.00 

198 

1.19e- 

4 

1.93 

258 

2.50e - 

1 

-0.00 

198 

1.19e-4 

1.93 

514 

2.50e- 1 

-0.00 

390 

3.04e - 

5 

1.97 

514 

2.50e - 

1 

-0.00 

390 

3.04e - 5 

1.97 

1026 

2.41e- 1 

0.05 

774 

7.69e- 

6 

1.98 

1026 

2.50e - 

1 

-0.00 

774 

7.70e - 6 

1.98 

2050 

2.30e - 1 

0.07 

1542 

1.93e- 

6 

1.99 

2050 

2.50e - 

1 

0.00 

1542 

1.94e- 6 

1.99 
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Figure 1: The relative errors in the norm versus the order of the matrix. The numbers in the legend 
are the values of the wavenumber. 


results of this experiment in Table They are also plotted in Figure These numerical results show 
that the OPGM are much more accurate than the CGM for all chosen k values. It is also illustrated in 
Fig. [T} In particular, when k = 5000 or 50000, the CGM totally fails to solve the equation while the 
OPGM still preserves the approximation accuracy. The convergence order of the OPGM is shown to be 
0{h?) in Table (c) and (d) when k is large enough which verifies the result of Theorem 

In experiment two, we solve equation with fixing the order of the coefficient matrices and varying 
the wavenumber from 10 to 10^. The orders of the coefficient matrices for the CGM and the OPGM are 
set to be 258 and 198, respectively. We plot in Fig the errors of both the methods with respect to the 
wavenumber and in Figj^the condition number of the coefficient matrix with respect to the wavenumber. 
From Fig we have two observations: The CGM is greatly affected by the wavenumber before the 
method fails to get an approximation with an acceptable accuracy; the OPGM nearly keep the uniform 
error with respect to the wavenumber when k is more than about 200. It verifies that the convergence 


6.2 
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Figure 2: The relative errors in the norm versus the wavenumber. The orders of the matrices are 258 
and 198 for the CGM and the OPGM, respectively. 
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Figure 3: Gondition numbers versus the wavenumber. The orders of the matrices are 258 and 198 for the 
GGM and the OPGM, respectively. 


of the OPGM is independent of the wavenumber. In Fig when k is larger than 300, the condition 
number of the coefficient matrix of the OPGM is uniformly bounded. This is consistent with that of the 
GGM. This conhrms the theoretical result in Theorem |7.2[ Though the GGM has the uniform condition 
number bound independent of k, it cannot give an appropriate approximation when k becomes large. 
This is because the standard spline spaces cannot approximate well an oscillatory function. In summary, 
the OPGM gives the optimal convergence order and at the same time preserves the uniform condition 
number of the GGM. 

To close this section, we note that to solve equation (2.1), we need to evaluate three kinds of integrals 
to obtain the linear system: (i) P £ Q (h) (/>P £ (hi) 

, Bj^ 2 e'^i^'^^), P £ j Q £ ^4 , where Bj ^2 is a B-spline basis function and Tp := p — 2, p e Z^. 
Integrals of class (i) can be evaluated explicitly. The other two classes may be computed numerically. In 
this paper, since K{s,t) = 1, {s,t) G P and the expression of / is simple, we can calculate the elements 
of the linear system through explicit expressions with the help of (4.15). However, for general K and /, 
we have to compute oscillatory integrals of classes (ii) and (hi) numerically. Highly oscillatory integrals 
have been understood deeply and calculated extremely efficiently by many methods, such as asymptotic 
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method, Filon-type method isniiQ], Levin-type method [32] , steepest descent method [25] and Clenshaw- 
Curtis-Filon-type method naimiis]. A new method is proposed in [29j by using the idea of graded 
mesh to analyze oscillatory integrals which is proved to be even more efficient than most of the existing 
methods. 
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